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Abstract 

We introduce a prototype model for globally-coupled oscillators in which each element is given 
an oscillation frequency and a preferential oscillation direction (polarization), both randomly dis- 
tributed. We found two collective transitions: to phase synchronization and to polarization order- 
ing. Introducing a global-phase and a polarization order parameters, we show that the transition to 
global-phase synchrony is found when the coupling overcomes a critical value and that polarization 
order enhancement can not take place before global-phase synchrony. We develop a self-consistent 
theory to determine both order parameters in good agreement with numerical results. 
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In recent years a considerable interest has been devoted to the self-organization properties 



1]. The work of Winfree 



first 



exhibited by networks of coupled nonlinear oscillators 
showed that the study of self-sustained non-identical oscillators is a suitable framework to 
achieve insight on the synchronization processes in biological systems. Based on Winfree's 



approach, Kuramoto [3| proposed a treatable mode 



for synchronizing oscillators successfully 



exploited in many fields, from heart physiology |4( to superconducting junctions j^. The 
underlying idea behind this success is that in many instances the dynamics of the individual 
oscillators can effectively by described as a limit cycle in which only one phase plays a 
relevant role. Then, for small disorder and weak coupling the Kuramoto model provides an 
excellent description of the synchronization process. A limitation of this model is that it 
does not consider the possible different direction of oscillation of the coupled oscillators. In 
fact the relationship between phase synchronization and a possible collective ordering of the 
oscillation direction has not been yet addressed. This question is of direct relevance in the 
field of optics: the cooperative behavior encountered in laser arrays has been investigated 
both from experimental [f], Q and theoretical jf], 0] points of view including descriptions in 
terms of the Kuramoto model [^] where the global coupling arises from light feedback from 
an external mirror. However, the vectorial nature of the electric field imposes a fundamental 
limitation to the description in terms of single phase oscillators. This description can only be 
used when the polarization degree of freedom is completely fixed by natural constrains. This 
is not the case, for example, in arrays of vertical-cavity surface-emitting lasers (VCSELs) 1(3 ] . 
where the polarization of the emitted light is not fixed by the structure ^j] , and the interplay 
between polarization and electric field global-phase requires at least a description in terms 
of two phases for each element. Indeed, it is possible to have states in which the global 
phases are synchronized despite of a misaligned polarization configuration. Such states 
have been observed experimentally in VCSEL arrays 13]. Moreover, polarization dynamics 
play an important role in the synchronization of master-slave VCSEL configurations, and 
polarization encoding has been recently proposed for high bit-rate encryption in optical 
communications jl^j ]. 

In this Letter, we develop an extension of the Kuramoto model as a prototype for the 
study of the fundamental properties of coupled oscillators described by vector fields in which 
at least two phases play a critical role: One associated with the natural oscillation frequency 
as in the Kuramoto model, and the other with the direction of oscillation (polarization). 
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We study the synchronization properties of an ensemble of globally coupled non-identical 
oscillators and show the existence of two transitions: phase synchronization and polarization 
direction ordering. We develop a self-consistent theory to determine the thresholds for both 
transitions and show that polarization ordering can never take place if the system is not 
already synchronized in frequency. 

Our analysis is made in the context of a general model, the Vector Complex Ginzburg- 
Landau Equation (VCGLE), which has been used for modellin g di fferent physical systems, 
from two-components Bose condensates to non-linear optics (l^ including laser emission 



12, UJ- The VCGLE can be written on 



from wide aperture resonators such as VCSELs 
symmetry grounds, but the determination of the parameters in the equation requires a 
specific physical model. We consider here parameter ranges of interest in optics. A set of N 
globally-coupled space-independent VCGLEs is given by 

d t Af = ( fij+iU j 3 )Af-(l + zf3 J )(\Af\ 2 + 1 ,\AJ\ 2 )Af 

N 
k=l 

where Aj (Aj) is the circularly right (left) polarized component of the j t h (j = 1 • • • N) vec- 
tor variable, uij is the natural oscillation frequency, (3j gives a nonlinear frequency shift and jj 
(a real number for lasers) couples the polarization components. The term (7 a +^7 P ) exp(±iSj) 
represents an external forcing [l7[ that linearly couples At and Aj. For example, for a 
VCSEL, the forcing arises from device anisotropies (dichroism and birefringence) [ll| that 
couple the circularly polarized components of the electric field, introducing two preferential 
polarization directions. Another example is a ring laser where any localized change in the 
refraction index breaks the invariance along the ring, introducing the same coupling be- 
tween the two counter-propagating modes and setting a preferential phase relationship 
between them. C is the strength of the global coupling which in laser arrays may be induced 
by external reflections (e.g. by placing a reflection at the common focus of the array j3|) or 
by a common active medium 13|. We introduce A^j = Q^j exp (itpf)- We consider jj < 1 . fo r 
which linearly polarized states (Q^ = Qj) are stable solutions of the solitary oscillators [12], 
as it is the case of VCSELs Jjj . Close to these solutions, we neglect the dynamics for each 
polarization component amplitude (Qf = 0), so that the system (JTJ can be described in 
terms of phase equations for each oscillator: The global phase <pj = (y?+ + ip~)/2, and the 
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rotational phase ipj = (<^~ — tpj )/2, whereas the latter determines the linear polarization 
direction. We have 

N 

X 

k=l 

\ 

A" 



C 

ipj = 7 a sin(2^j - Sj) + — ^2 sin (V'fe - i>j) cos {<t>k ~ <f>j) , (2) 

k=l 



c 

4>j = LOj + 7 P cos(2tAj - Sj) + — 22 sin (4>k - 4>j) cos (ip k - Ipj) . (3) 

fc=l 



In the uncoupled case (C = 0), the global phases 0j rotate at a constant frequency, 
whereas the polarization angles ipj reach a steady state, thus modelling a solitary laser 
emission. In fact, for C = we have two orthogonal linearly polarized solutions for the jth 
oscillator: 2ipj = Sj, <pj = <poj + (ujj + J p )t, and 2ipj = Sj + n, <pj = <p j + (ujj — J p )t, where 
cpoj is a constant. For 7 a < the first solution is selected, whereas for 7 a > the second 
is selected. In laser physics, the parameter 7 a models the different linear gain encountered 
by the two linearly polarized solutions, thus making linearly stable the solution with the 
ligher linear gain. In the same context, the parameter 7 p models the cavity birefringence 



which splits the emission frequency of the two orthogonal linearly polarized solutions 
by an amount equal to 27 P . In the following, we take 7 a < 0, so we will refer to 2ipj = Sj 
as to the natural polarization angle of each oscillator. Our results, however, do not depend 
on this choice, neither do on the sign of j p , which we set positive. Fixing the polarization 
degree of freedom (2ipj(t) = Sj = S for all j) equation (jSJ) reverts to the Kuramoto model: 

<t>3 = u i + 1p + § £*Li sin (^ _ 

The differences in the natural polarization angles and frequencies of the oscillators repre- 
sent two sources of disorder in our system. They are statistical quantities, randomly chosen 
from two symmetric unimodal distributions q(S) and p(u), with zero mean and standard de- 
viation a u and as, respectively. Therefore, we introduce two order parameters to characterize 
the degree of phase synchronization and polarization ordering, respectively 

1 - 

V eX P (*X) = ^ ' 

k=l 
N 

X 

k=l 

Without coupling p averages to zero while, as 2ipj — Sj, rj accounts for the natural disorder 



p exp (id) = — Y ex P (^*) • ( 5 ) 



in the polarization angle. In the continuum limit, rj = r] 



n n ^p(io)q(S)dS 



which is 



non-zero unless q(S) is a uniform distribution between —it and 7r. 



For small coupling the global phases (f)j are de-synchronized, which leads the coupling 
term in the polarization Eq.Qto the vanish. Therefore each oscillator remains oscillating 
in its natural polarization angle. No polarization interaction takes place until the phases 4>j 
start to synchronize. Increasing C, two different scenarios toward polarization ordering and 
phase synchrony (i]=p=l) are found depending on the relative strength of the polarization 
as and phase disorder cr w . 

For o~ w C(7i, the transitions to phase and polarization synchrony are well separated. The 
phases <pj synchronize first. The transition to phase synchrony can be analyzed by taking 
2ipj = Sj (frozen polarizations), so that the set 0-© can be approximated by 

h = Uj + lp + jz sin (<Pk ~ <f>j) cos ( k 3 J . (6) 
k=i \ ' 

Averaging the polarization angles, reduces Eq. © to a Kuramoto-like model with an effective 
coupling C 

C N 

<fij = Uj + % + Jr Sil1 (0* - <t>j) > ( 7 ) 

k=l 

where C = C J cos [(5 — S')/2] q(S)q(S')dS = Crfo. The polarization disorder makes the phase 
coupling less effective but not vanishing. Following the standard treatment of the Kuramoto 
model [3|, the self-consistent equation for the order parameter amplitude p reads 

p = Cp cos 2 {4>)q{Cpsm((f)))d(p. (8) 

J-n/2 

Therefore, the critical coupling C t for the onset of collective phase synchronization reads 

2 

* = np(0) J cos [(5 - 6')/2] q(6)q(6')dS ' (9) 
Fig. ^ shows the good agreement between the transition to phase synchronization obtained 
from numerical integration of ©-© and the solution of the self-consistent Eq. (JHJ). Notice 
the excellent agreement obtained for the onset of synchronization given by Q, Ct = 0.01968. 
The distribution of averaged dressed frequencies Q = (^<pj (left inset of Fig. [J) shows a highly 
dominant peak which comes from the synchronized oscillators (0.97iV in this case). Notice 
also that for C < Ct the polarization order parameter takes a constant value rjo which 
corresponds to the initial polarization disorder, in agreement with the assumptions leading 
to ©. 



Increasing further the coupling strength, the oscillators leave the respective natural polar- 
ization angles and start to order in polarization. As the phase synchronization has already 
been achieved, we are now in the position to develop a self-consistent theory for the po- 
larization ordering as follows: Assuming perfect phase synchronization = <j)j), Eq. (J2J) 
becomes 

C N 

fa = 7a sin(2^ - 8j) + — sin (Vfc - *l>j) ■ ( 10 ) 



N 
fc=i 



Since the individual polarization is not a self-oscillating dynamics, Eq. (|10|) is not a 



Kuramoto-like model. However, from (0J) we have r]sm(x — = h X)fe=i srn (^fc — V^X 
which introduced in Eq. (jl(Jj) yields 

ipj = -f a sm(2^- - Sj) - Cr) sin('0j - x) ■ (H) 

The stationary solution i>j(6,r],x), given implicitly by 

7 a sin(2^ - 5,-) - CV/ sin(^ - x) = , (12) 

can be introduced in Eq. (jlj) to self-consistently find r] and x- I n the continuum limit we 
have 

r]exp(ix) = J exp(i$(5,r},x))q{S)d5 . (13) 

Altogether, Eqs. (|T2"j) and ()13|) allow for the calculation of the polarization order parameter, 
for example through a Newton-Raphson method, so that the polarization ordering can be 
fully described. The imaginary part of integral (|13|) was found to vanish (x = 0), if q(5) 
is even. Fig. shows the agreement between the evaluation of r] using the definition (JJJ) 
with the results of the numerical integration of the full set Eqs. (|2*|-(|3*| and using the self- 
consistent approximation given by Eqs. ([T2^ - ([TH)l . We obtain good results even for small 
coupling where global phases are desynchronized (in that regime the contribution of the 
coupling term in the polarization equation is negligible). The ordering of the polarization 
directions induces a loss of coherence where the phases partially de-synchronize, lowering 
p. The reason is that, as the polarization order is increased, the polarization angles depart 
from the natural angle, and therefore the term 7 P cos (2^- — Sj) in Eq. (JHJ) plays the role of an 
added disorder to the natural frequencies uij. Increasing j p , this effect is linearly increased, 
enhancing the coherence loss extent, as shown in Fig. ^ For 7 p = 5, p is reduced down 
to 0.65. The averaged dressed frequency distribution (right inset of Fig. |TJ) shows that the 



peak at = is lowered in the same proportion and two lateral lobes associated to drifting 
oscillators appear, yielding an overall shape for the distribution similar to that of partially 
synchronized Kuramoto oscillators [3] . From a practical point of view, the coherence lowering 
would have a direct impact the output intensity in VCSEL arrays. A reduction of p down 
to 0.65 leads to coherent output intensity of only 40% with respect to the fully synchronized 
case. Finally, for large coupling, complete phase synchronization and polarization ordering 
are achieved. 

Numerical simulations for different values of as showed that decreasing the disorder in 
the natural polarization angles, the polarization transition to synchronization shifts to lower 
values of the coupling. However, a polarization order enhancement is not possible before 
the phases start to synchronize, so for a u ~ as or a^ > as the two transitions take place 
simultaneously. Nevertheless, the polarizations are still effectively uncoupled until the phases 
start to synchronize, so the self consistent equation (jHJ) still holds as well as the prediction 
for the phase synchronization onset Ct, which now also signals the onset of the polarization 
ordering as shown in Fig. [21 The self-consistent equation for the polarization order parameter 
still gives a good description of the polarization order enhancement. 

In conclusion, we have introduced a theoretical framework to study the synchronization 
properties of a system of globally coupled oscillators extending the results for limit cycle 
oscillators to include the consideration of oscillation direction (polarization). Two sources 
of disorder are included: Randomly distributed natural frequencies and natural oscillation 
directions. Increasing the coupling no polarization order enhancement is possible until the 
phases start to synchronize, because the phase disorder destroys the interaction among the 
polarization variables. This is in agreement with experimental results observed in VCSEL 
arrays Q|. Typically, the frequencies synchronize first, and polarization synchrony takes 
place at a higher coupling level, through a partial de-synchronization of the phases (coherence 
lowering). We have developed self-consistent approximations which provide a very good 
estimation of the synchronization properties of system. Increasing the disorder in the natural 
frequencies or decreasing the disorder in the natural polarization angle the two transitions 
merge in a unique process to full synchrony, and we provided the critical coupling for its 
onset. 
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FIG. 1: Amplitude of order parameters p and rj as function of the coupling C . The solid line 
corresponds to the theoretical predictions, whereas the symbols to numerical integration of © and 
(j3J) for typical VCSEL birefringence values (o: 7 p = 0, +: j p = 1, x: 7 P = 2.5, A: j p = 5). We 
have considered 7 a = —0.5, ./V = 10 3 , a Gaussian distribution p(u>) for the natural frequencies with 
Gu = 10~ 2 , and a uniform distribution q{5) = ^ for —A < 5 < A, with A = tt/2 (as = = 
0.9068) for the natural polarization angles. The insets show the time-averaged dressed frequencies 
distribution for 7 p = 5, C = 0.1 (left) and C = 2.5 (right). The natural frequency distribution 
pipS) is shown for reference (solid line). 
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FIG. 2: Amplitude of order parameters p and r\ as in Fig. but with larger natural frequency 
disorder a u = 0.12, (now C t = 0.2362). 
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